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Abstract 



The electromigration (EM) of metallic lines is studied in terms of competition be- 
tween two percolative processes taking place in a random resistor network. The 
I effects associated with the transport of mass and with the consequent growth of 

' the internal stress are accounted for by stochastic generation and recovery of voids, 

\^ . driven by the external current. Monte Carlo simulations enable us to investigate 

within a unified theoretical framework a variety of relevant aspects of EM degra- 
dation as: damage patterns, the distribution of the time to failures, early stage 
. resistance change in the presence of compositional effects. Black's law, geometrical 

^ , effects, etc. A general agreement is found with the experimental findings. 

a 

a 
o 
o 

> 

■ In recent times, the study of breakdown phenomena in electrical devices has 

received a lot of attention, mainly for the increasing level of miniaturization 
which makes failure problems more and more crucial [1]. Failure occurs in 
many cases by the degradation of metallic interconnects (thin films) which, 
because of electromigration (EM), lose their conducting properties. EM arises 
from the transport of matter at the atomic level driven by a high current den- 
sity [1-3] . More precisely, the damage, consisting in the formation and growth 
of voids and hillocks in different regions of the film, is due to a non-steady 
atomic trasport associated with a nonvanishing divergence of the atomic flux. 
A fundamental ingredient in the description of the EM damage of the inter- 
connects in electronic devices is represented by the granular structure of the 
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materials employed, Al, Cu, Al alloys, etc. In fact, it has been observed [1-3] 
that the atomic transport through the grain boundaries (transport channels) 
far exceeds that through the grain bulks. Therefore, it is generally possible 
to neglect mass transport everywhere except within these channels and to de- 
scribe the film as an interconnected grain boundary network [1]. Puthermore, 
a high degree of disorder is usually present in alloy films due to compositional 
effects (CE) and to thermal gradients [1-3]. Finally, it must be underlined 
the role of mechanical stress. Indeed, the depletion and accumulation of mass 
in different regions of the film, implies also the growth of mechanical stress, 
which contrasts the atomic flux due to EM [1,2,4]. On the basis of the above 
description, we have developed a model which describes the EM damage of 
interconnects in terms of two competing biased percolations taking place in 
a Random Resistor Network (RRN) [5]. Monte Carlo simulations allow us to 
explore a variety of relevant aspects of EM degradation as: damage patterns, 
the distribution of the time to failures, early stage resistance change in the 
presence of compositional effects. Black's law, geometrical effects, etc. The 
paper is organized as follows. Section 2 describes the model used and the nu- 
merical procedure of simulations. Section 3 reports the results and compares 
theory with well accepted phenomenological behaviors. Major conclusions are 
drawn in Sec. 4. 



2 Model 



A metaUic fine of length L and width w is described as a two-dimensional 
rectangular-lattice network of regular resistors, r^eg- The network lies on an 
insulating substrate at temperature Tq acting as a thermal reservoir. The total 
number of elementary resistors in the network is given by: Ntot = '^NlN^ + 
Nl—Nw, where Nl and determine the length and the width of the network. 
To reduce the computational efforts, long lines are simulated by assuming 
that only a portion of the line is responsible of the resistance variations. The 
network represents this portion and its length is taken 1/F times smaller 
than the line length. Thus, the relative resistance variation of the whole line 
is obtained by multiplying the relative resistance variations of the network 
by the factor F. The resistance of the n-th resistor depends on temperature 
according with: 

TregATn) = Trefi^ + Oi{Tn - T^e/)] (1) 



where a is the temperature coefficient of resistance (TCR), T„ is the local 
temperature, T^e/ and r^e/ are the reference values for the TCR. The network is 
contacted at the left and right hand sides to perfectly conducting bars through 
which a constant stress current / is applied. The Joule heating induced by / 
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is taken into account by defining T„ as [6]: 



J-^neig rn=l 



where A is the thermal resistance of the single resistor and Nneig the number 
of first neighbours of the n-th resistor. The value S = 3/4 is chosen to provide 
a uniform heating for the perfect network [6] . By taking this expression we are 
assuming an istantancous thcrmalization. In the presence of a current stress, 
the EM damage is simulated by allowing the breaking of regular resistors, 
i.e. by replacing Vreg ^opi ^op = 10^ r^e/- We remark that rop represents 
an open circuit (OP) associated with the formation of voids inside the line. 
Simultaneously, a recovery process which allows the healing of the broken 
resistors, is also considered. The complete failure of the line is thus associated 
with the existence of at least one continuous path of OP between the upper 
and lower sides of the network (percolation threshold) [5]. The probability 
Wop for the n-th resistor to become OP is taken as: 

Wop = eM-Eop/ksTn) (3) 



where Eop is a characteristic activation energy and ks the Boltzmann con- 
stant. A similar expression is taken for the recovery process probability Wre, 
with Eop replaced by Ere- Accordingly, the elimination of resistors represents 
the mass transport due to an electron wind, and the recovery of resistors, the 
back-up fiow of atoms due to the mechanical stress. The two activation ener- 
gies control the competition between these two opposite processes. 



2.1 Compositional effects (CE) 



In Al alloys (Al-Cu, Al-Si) CE are observed at the early stages of EM test 
of metallic lines. Here, we limit ourself to study the CE due to: (i) a sample 
heating, (ii) an intrinsic defectiveness of the manifactured samples. To sim- 
ulate the CE induced by heating, we take into account the actual thermal 
treatment undergone by the lines just before the EM test. Here, we consider 
one specific EM experiment which is performed with standard techniques on 
Al-0.5%Cu lines. For these lines, the last thermal treatment occurred during 
fabrication was a high temperature annealing followed by a rapid cooling. This 
treatment left a non-equilibrium concentration of Cu dissolved into the Al ma- 
trix. Therefore, in the early stage of the EM test, the heating associated with 
the stress conditions allows a relaxation toward the equilibrium concentration 
of the Cu atoms and gives rise to the formation of CuAl2 precipitates. The 
variation of the line composition is simulated by changing, with probability 
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Fig. 1. Flowchart of the Monte Carlo simulator. 

Wri, the resistance of the n-th resistor from rreg to rimp (impurity resistors), 
where we take: 

Wri = exp{-Eri/kBTn) (4) 

being Ej-i an activation energy characterizing this process. The mechanism of 
Cu dissolution into the Al matrix, antagonist to the precipitation, is simul- 
taneously considered by taking a similar expression for the probability of the 
reverse process Vimp — >■ Vreg, which is determined by an activation energy Eir. 
As a consequence of the two competing processes, a steady state concentration 
of Vimp is achieved on a time scale much shorter than that of EM. The intrinsic 
defectiveness of the samples manifests itself by the existence of samples of the 
same material and geometry but with drastically different failure times. To 
simulate this second case of CE, we have considered the evolution of networks 
with different initial concentration pini of OP resistors. Accordingly, we have 
assumed pini randomly distributed, and chosen values of ro in such a way to 
keep nearly the same resistance value of similar undefected samples. 

2.2 Numerical procedure 

The network evolution is obtained by Monte Carlo simulations which are car- 
ried out according with the flowchart shown in Fig. 1. (i) Starting from the 
initial network, we calculate in and the network resistance R by solving Kirch- 
hoff's loop equations. Moreover, we calculate T„ by using Eq. (2). (ii) OP and 
Timp are generated with the corresponding probabilities Wqp and Wri, while 
the remaining Vreg are changed according to T„. Then, in and T„ are recalcu- 
lated, (iii) OP and r^mp are recovered, (iv) i„, T„ and it! are recalculated. This 
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Fig. 2. Simulated damage pattern at the complete failure of a network 40 x 40 by 
using the parameters specified in the text. 

procedure is iterated from (ii), thus the loop (ii)-(iv) corresponds to an itera- 
tion step, which is associated with a unit time step on an arbitrary time scale 
to be calibrated by comparison with experiments. The iteration proceeds until, 
depending on the values of the model parameters, the following two possibili- 
ties are achieved: irreversible failure or steady-state evolution. The irreversible 
failure can be established by a criterion for the increase of resistance over its 
initial value as suggested by experiments (typically of 20 %). 



3 Results 

To check the ability of the model to describe EM experiments, simulations are 
performed on Al-0.5%Cu lines [7]. A significant aspect related to the electrical 
failure is the damage pattern, which evidences the correlations among different 
failing regions just before the complete failure. Figure 2 reports a typical 
damaged pattern obtained for a network of sizes 40 x 40 with the following 
parameters chosen according to the experimental conditions of Ref. [7]: / = 
1.1 X 10-2 A, To = 219 °C, Tref = °C, rref = 296 n, a = 3.6 x 10"=^ "C-\ 
To provide for an initial heating of the network of about 15 °C, the value of 
the thermal resistance A is taken as: 1.5 x 10^ °C/W; to account for an initial 
defectiveness of the line we take pini ~ 0.01. Moreover, to save computational 
time, we choose F=30 and for the activation energies responsible of EM we use 
values smaller than those found in experiments, as: Eqp = EptE = 0.43 eV. 
We have tested that, by scaling the value of Epp in the range 0.2^0.6 eV, the 
median time to failure (MTF) changes over more than two orders of magnitude 
while its relative standard deviation remains constant within a factor of two. 
Consistently with this choice, to account for the different characteristic times 
of the two processes, also the activation energies responsible of CE have been 
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Fig. 3. Simulation of the relative variations of resistance as function of time for 
AI-0.5%Cu lines. The insert shows enlarged CE in the early stage of EM. 



reduced: we take Eri — 0.22 eV and £^^^==0.17 eV. Moreover, for the impurity 
resistor we use rimp — 350 Jl. We note that the damage pattern in Fig. 2 well 
reproduces the experimental pattern observed by SEM in metallic lines [1,3]. 
The results of simulation for the resistance evolutions arc reported in Fig.3. 
All the parameters are the same used for Fig. 2, except for the network sizes, 
that are now Nl — — 70. The insert shows the early stage resistance 
decrease due to the thermal CE. The simulations reproduce satisfactorily the 
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Fig. 4. Simulation of the relative variations of resistance as function of time corre- 
sponding to networks with different initial defectiveness and comparable resistance. 
Continuous curves represent evolutions obtained by taking Pini « 0.1, the dotted 
curves are for pini « 0.3. The insert reports a set of experimental results performed 
on samples with the same geometry and comparable resistance. 



6 



experimental evidence of both the early stage and the pre-failure regime of the 
resistance evolution [7] . To investigate the role of the intrinsic defectiveness of 
the line, we have simulated the evolution of the network with a different initial 
concentration of OP resistors randomly distributed. Moreover, we have chosen 
rref — 229 Q to have nearly the same resistance of the networks considered in 
Fig. 3. We report in Fig. 4 the relative variations of the resistance as function 
of time for the two sets of simulations: continuous curves are the same of 
Fig. 3 {pini ~ 0.1), the other ones corresponding to Pini ^ 0.3. We show in 
the insert the results of the experiments [7] performed on samples with the 
same geometry and comparable resistance but exhibiting completely different 
failure times. The simulations well reproduce the experiments [7]. 

EM tests are performed in the so-called "accelerate conditions" , i.e. the cur- 
rent density, J, and the stress temperature, T, are much higher than those 
corresponding to normal operation conditions. This is an obliged choice to 
obtain results on a time scale of days or weeks, whereas, in normal conditions, 
the failure times are of the order of years. Therefore, the evaluation of the fail- 
ure times in normal conditions, comes from the existence of a scaling relation 
connecting accelerate conditions with normal ones. This relation is provided 
by the well known Black's law [8] for the median time to failure (MTF): 



MTF-^ oc J" T-^ exp 



-E 



where E is the activation energy for EM. For moderate stress conditions the 
value of 11 is found between 1 and 2, while higher values are found for extreme 
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Fig. 5. Median time to failure as a function of stressing current for a rectangular 
network 12 x 400 with the parameter values specified in the text. The slopes in the 
relevant regions of the curve (low, moderate, and high current stresses) are reported 
in the figure. 
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stress conditions. Several models have been proposed to explain this law with 
contradictory predictions concerning the value of the exponent n [9-11]. There- 
fore, we have calculated the MTF as a function of the applied current and the 
results are shown in Fig. 5. The simulations are performed on rectangular net- 
works with sizes 12 x 400 by taking the same values of To, T^e/, « used before 
and the results are shown in Fig. 5. The values of the other parameters are 
the following: r^e/ = 0.048 Vt, n^p = 0.016 Vl, A = 2.7 x 10^ "C/W, pini = 0, 
F = 200, Eop = 0.41 eV, Ere = 0.35 eV, = 0.22 eV and = 0.17 
eV. From Fig. 5 we can see that experimental behaviors are naturally repro- 
duced by our simulations, i.e. n ~ 1 for low bias, and n grows up to 5.8 for 
extreme bias values. It is interesting to consider the distribution of the times 
to failure (TTF) of a set of identical samples in the same stress conditions. 
Therefore we report in Figure 6 the distribution of the TTF obtained by sim- 
ulations, at increasing values of the stressing current (going from the right to 
the left). All the parameters are the same of that used for Fig. 5. For each 
current value, 20 realizations have been considered. Precisely, in the figure the 
cumulative failure percentage is reported as a function of the logarithm of the 
TTF and a transformation has been applied to the ordinate axis in such a 
manner that a log-normal distribution would be a straight line. From the fig- 
ure we can conclude that the distribution of the simulated TTF is log-normal 
(within errors) and that the shape factor (the log-normal standard deviation) 
is independent of the applied current. These conclusions well agree with the 
experiments [1-3]. As a result of the competition between the electromigration 
and the backfiow of mass determined by stress gradients, the failure can occur 
or not depending on the length of the line. As a consequence, it is possible to 
consider this geometrical effects on the value of the MTF. Figure 7 shows the 
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Fig. 6. Distribution of TTFs obtained from 20 simulations of failure at different 
stressing currents. The network sizes and the parameter values are the same of Fig. 
5. Lines show the fit with a log-normal distribution. 
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Fig. 7. Median time to failure as a function of the length to width ratio at fixed 
width. The network width and the parameter values are the same of Fig. 5. Points 
are the results of average over 20 failure simulations. 

simulated MTF as function of the length to width ratio of the network (with 
fixed width). We remark that the MTF diverges for small lengths, according 
with the experimental results [2]. Experiments also show that for a line of 
given length it exists a critical value of the current density below which there 
is no electromigration. Moreover, this breakdown value is inversely propor- 
tional to the length [1-3]. Our simulations well reproduce also this behavior 
as shown in Fig. 8, where the breakdown current is reported as a function of 
the width to length ratio (again at width fixed). 
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Fig. 8. Breakdown current as a function of the width to length ratio at fixed width. 
The network width and the parameter values are the same of Fig. 5. The curve is a 
linear fit to data. 
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4 Conclusions 



We have developed a Monte Carlo approach to electrical breakdown phenom- 
ena, with particular attention on its apphcations to the study of electromigra- 
tion in metalhc hues. Simulations are performed on a RRN whose configura- 
tion is evolving because of two competing percolations driven by an external 
current. The results of the simulations well agree with the experimental find- 
ings concerning several features of the electromigration phenomenon, including 
compositional and geometrical effects. The flexibility of the approach allows 
wide possibilities of further extcntions of the model. Future goals are the 
extraction of information about precursor effects of the electical breakdown 
phenomena. 
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